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Using  empirical  models,  parameters  have  to  be  estimated  from  experimental  data.  Experimental  charac¬ 
terization  of  fuel  cell  stacks  is  an  expensive  and  time-consuming  task.  Therefore  it  is  very  important  to 
choose  an  experimental  design,  which  maximizes  the  statistical  quality  of  the  resulting  information.  Box 
and  Lucas  (Biometrika  46  (1959))  showed  that  it  is  possible  to  optimize  nonlinear  experimental  designs 
by  the  minimization  of  the  covariance  matrix  of  the  least  squares  estimate.  The  aim  of  this  work  is  to  adopt 
this  general  method  in  order  to  investigate  its  ability  for  application  in  polymer-electrolyte-membrane 
fuel  cell  (PEMFC)  characterization.  Based  on  an  empirical  PEMFC  model  a  D-optimal  design  criterion 
has  been  developed  and  validated.  Numerical  methods,  evolutionary  and  heuristic  are  investigated  with 
respect  to  fast  and  robust  evaluation  of  the  design  criterion.  For  a  given  set  of  experimental  data  best 
results  are  achieved  using  a  heuristic  approach,  a  so-called  sequential  search.  Based  on  that  result 
an  algorithm  to  obtain  an  optimal  design  of  experiments  (DoE)  in  a  nondeterministic  operating  area 
is  introduced.  The  proposed  algorithm  is  able  to  take  into  account  experimental  limitations  due  to 
test  facilities  or  examinee.  The  algorithm  further  allows  to  include  existing  and  for  reference  needed 
experiments. 

©  2008  Elsevier  B.V.  All  rights  reserved. 


1.  Introduction 

Intelligent  energy  management  is  a  cost  effective  key  path  to 
realize  efficient  automotive  drive  trains.  To  develop  operating  strat¬ 
egy  in  fuel  cell  drive  trains,  precise  and  computational  efficient 
models  of  all  system  components,  especially  the  fuel  cell  stack,  are 
needed. 

System  identification  is  an  essential  step  in  empirical  modelling. 
Especially,  the  estimation  of  unknown  parameters  is  a  typical  prob¬ 
lem  in  the  development  of  a  PEMFC  model.  The  standard  method 
is  the  analysis  of  experimental  data  from  measurements.  Due  to 
system  complexity  of  a  PEMFC  the  experimental  investigation  is 
an  expensive  and  time-consuming  task.  In  most  of  the  cases  just 
a  limited  operating  area  is  available  for  the  measurements  caused 
by  the  fuel  cell  capabilities,  physical  laws,  test  equipment,  or  inter¬ 
actions  between  these  components.  To  give  a  simple  example  for 
equipment  limitations  just  consider  the  gas  supply  to  the  cell.  Mass 
flow  in  general  is  limited  to  a  certain  range  and  therefore  the  sto¬ 
ichiometry  X  {X  =  ratio  of  supplied  reactant  to  chemically  needed 
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reactant)  is  also  limited  in  dependency  of  the  current  density.  At 
a  low  load  point  using  maximum  gas  flow  a  higher  stoichiometry 
can  be  realized  than  at  high  load  points.  Also  the  fuel  cell  itself  will 
work  properly  with  low  gas  pressure  at  low  loads,  but  not  at  high 
loads. 

Therefore  it  is  on  the  one  hand  important  to  minimize  the  exper¬ 
imental  work  with  consideration  of  these  constraints  and  on  the 
other  hand  to  achieve  the  best  possible  parameter  estimation.  The 
aim  of  this  article  is  to  show  and  demonstrate  a  systematic  method 
for  investigation  of  a  fuel  cell,  which  performs  both  tasks  simulta¬ 
neously.  It  is  known  that  the  quality  of  the  parameter  estimation 
depends  on  the  used  measurement  points  and  hence  on  the  design 
of  experiments  (DoE).  Therefore  the  issue  can  also  be  understood  as 
to  find  the  data  set,  which  delivers  the  largest  information  content 
if  experimental  data  are  available  or  if  not  to  provide  an  optimized 
experimental  design  for  a  PEMFC  stack. 

Box  and  Lucas  [1  ]  show  in  their  work  a  mathematical  theory  for 
increasing  the  quality  of  least  squares  based  parameter  estimation 
by  minimizing  the  covariance  matrix  of  its  solution.  Several  criteria 
are  published  in  the  literature  [2]  to  evaluate  this  optimization.  In 
this  paper  these  criteria  are  compared  and  the  most  promising  one 
is  applied  on  a  set  of  405  test  points,  to  find  the  set  of  50  points 
which  is  most  suitable  for  the  parameter  estimation  of  an  empirical 
fuel  cell  voltage  model.  Moreover,  it  is  investigated  whether  this 
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Table  1 

Typical  factor  levels  of  selected  inputs  in  PEMFC  characterization 


Table  2 

Experimental  conditions 


State 

Symbol 

Operating  range 

Cathode 

Anode 

Pressure 

Pc,  Pa 

1.05-2.5  bar 

1.05-2.5  bar 

Relative  humidity 

rhc,  rha 

0.32-1 

0.32-1 

Temperature 

T 

75  °C 

Current  density 

i 

0-2  A/cm2 

optimized  data  set  increases  the  quality  of  parameter  estimation 
and  therefore  the  stack  model  accuracy. 

2.  Design  of  experiments 

2.1  Design  methods 

If  the  operating  area  is  unrestricted  and  the  number  of  fuel  cell 
inputs  are  low,  a  full  factorial  DoE  can  be  performed.  Interactions 
as  well  as  nonlinearities  can  be  detected  if  more  than  three  factor 
levels  of  each  input  are  applied.  But  with  a  high  number  of  inputs  or 
factor  levels  the  experimental  effort  rises  exponentially.  In  Table  1 
the  typical  number  of  factor  levels  is  given.  A  full  factorial  design  of 
experiments  with  the  maximum  values  of  each  input  would  thus 
contain  25.920  experimental  points. 

If  a  system  behaves  linear  within  an  unrestricted  operation  area 
a  fractional  factorial  DoE  can  be  chosen.  That  technique  reduces  the 
number  of  experiments  compared  to  full  factorial  designs  but  is  not 
able  to  detect  interactions  between  inputs  or  nonlinear  behaviour. 
Only  a  fraction  of  the  edges  of  the  hypercube,  defined  by  the  dimen¬ 
sion  of  the  operating  area,  is  tested. 

A  trade-off  between  effort  and  ability  to  detect  interactions  and 
nonlinear  behaviour  gives  a  so-called  central  composed  design.  All 
edges  of  the  hypercube  are  used.  Additional,  centre  points  and  star 
points  are  added.  These  three  designs  are  easy  to  construct,  but  all 
require  an  unrestricted  operating  area  and  are  therefore  not  suited 
for  global  fuel  cell  characterization. 

A  first  approach  to  design  experiments  in  restricted  areas  is  the 
so-called  Box-Benken  design.  In  this  design  the  number  of  experi¬ 
ments  is  fixed  and  cannot  be  varied.  The  experiments  are  not  placed 
in  the  edges  of  the  operating  area.  More  details  about  the  men¬ 
tioned  and  additional  designs  can  be  found  in  DoE  literature,  like 
Refs.  [3-5].  More  suited  to  the  requirements  for  PEMFC  characteri¬ 
zation  are  optimal  DoEs. 


of  optimal  criteria  is  able  to  design  nonlinear  DoEs  with  a  low  num¬ 
ber  of  experiments  in  arbitrary  restricted  operating  areas  without 
being  limited  to  a  fixed  number  of  experiments  or  factor  levels. 

Ref.  [1  ]  shows  that  the  covariance  matrix  COV  of  the  estimated 
model  parameters  §  obtained  by  a  least  squares  estimation  can  be 
approximated  by  the  following  equation: 

COV  =  (FIM/  •  FIM)-1  a2.  (1) 

with  the,  in  general  unknown,  variance  a 2  and  FIM  the  Fisher  infor- 
mation  matrix.  The  FIM  is  a  <2>-by -0  matrix  of  the  partial  derivative 
fim  (Eq.  (2))  of  the  0th  model  parameter  for  the  dth  set  of 
experimental  conditions  ®d  (Eq.  (3)).  Here  d  is  the  run  index  of 
the  experimental  conditions.  It  starts  with  1  and  ends  at  &.  Addi¬ 
tionally,  6  is  the  index  of  the  model  parameters  in  the  range  from 
1  tO  0. 

FIM  =  {fim0d},  0  =  1  ...0,  d  =  1...0.  (2) 


filled 


9fim(6>,  £) 
9? e 


(3) 


The  basic  idea  of  all  optimal  DoEs  is  the  minimization  of  a  quality 
criterion  J  related  to  COV  or  (FIM'-FIM)-1  respectively.  This  opti¬ 
mization  can  be  achieved  by  several  design  criteria  [5]. 

The  trace-criterion  (A-criterion)  (Eq.  (4))  minimizes  the  trace, 


JA  =  trace 


(fimtfim)  1 


(4) 


the  determinant-criterion  (D-criterion)  (Eq.  (5))  minimizes  the 
determinant, 


Jd  =  det 


(fimtfim)' 


(5) 


2.2.  Optimal  DoE 

Since  powerful  computing  has  become  affordable,  more  com¬ 
plex  and  powerful  design  criteria  have  been  developed.  The  group 


Fig.  1.  Graphic  interpretation  of  the  optimal  DoE  criteria  for  minimizing  the  covari¬ 
ance  matrix. 


the  eigenvalue-criterion  (E-criterion)  (Eq.  (6))  minimizes  the  max¬ 
imum  eigenvalue  (EIG), 


JE  =  max  \  EIG 


|  EIG  (FIMtFIM)  1  ) 


(6) 


Table  3 

Experimental  details  of  campaigns  1-3  (rha/c  =  0.92) 


Campaign 

1 

2 

3 

Gas  pressure 

Pc  (bar) 

Pa  (bar) 

Pc  (bar) 

Pa  (bar) 

Pc  (bar) 

Pa  (bar) 

1.05 

1.05 

1.05 

1.50 

1.50 

1.05 

1.10 

1.10 

1.10 

1.50 

1.50 

1.10 

1.25 

1.25 

1.25 

1.50 

1.50 

1.25 

1.50 

1.50 

1.50 

1.50 

1.50 

1.50 

1.75 

1.75 

1.75 

1.50 

1.50 

1.75 

2.00 

2.00 

2.00 

1.50 

1.50 

2.00 

2.50 

2.50 

2.50 

1.50 

1.50 

2.50 

50 
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and  the  confidence-interval-criterion  (M-criterion)  (Eq.  (7))  mini¬ 
mizes  the  maximum  element  of  the  main  diagonal. 


JM  =  max  < 


(fimt-fim)~ 


Aoo 


0  =  1,2. ..</>. 


(7) 


Also  a  geometric  interpretation  of  these  criteria  is  possible,  whereas 
the  A-criterion  equals  the  mean  axis  length  and  the  E-criterion 
the  longest  axis  of  the  confidential  ellipsoid  spanned  by  the  esti¬ 
mated  parameters.  The  M-criterion  can  be  explained  by  the  longest 
side  of  an  imaginary  box  around  the  ellipsoid  and  the  D-criterion 
can  be  interpreted  as  the  volume  of  the  confidential  ellipsoid.  An 
illustration  of  these  criteria  is  drawn  in  Fig.  1. 


Table  4 

Experimental  details  of  campaigns  4-6  (pc/a  =  1.5  bar) 


Campaign 

4 

5 

6 

Gas  relative  humidity 

rhc 

rha 

rhc 

rha 

rhc 

rhc 

0.32 

0.32 

0.32 

0.92 

0.92 

0.32 

0.40 

0.40 

0.40 

0.92 

0.40 

0.80 

0.80 

0.80 

0.92 

0.80 

0.92 

0.92 

0.92 

0.92 

0.92 

1.00 

1.00 

1.00 

0.92 

0.92 

1.00 

3.  Experimental  data  and  fuel  cell  model 

3  A.  Experimental  data 

For  the  validation  of  the  investigated  design  criteria  a  set  of 
experiments  is  used.  In  total  405  different  data  points  are  available, 
since  33  experimental  runs  have  been  carried  out  varying  cathode 
pressure  pc,  anode  pressure  pa,  cathode  gas  relative  humidity  rhc 
and  anode  gas  relative  humidity  rha  and  a  complete  polarization 
curve  was  measured  for  each  run.  All  measurements  were  per¬ 
formed  under  isothermal  conditions  of  75  °C.  To  avoid  the  influence 


of  gas  supply  to  fuel  cell  voltage,  the  cell  was  fed  with  30  times  more 
reactants  than  the  cell  consumed.  An  overview  of  all  experimental 
conditions  is  given  in  Table  2. 

During  the  first  three  campaigns,  gas  pressures  are  varied  while 
relative  humidity  of  anode  and  cathode  gas  is  kept  constant  at 
rha/c  =  0.92.  During  campaigns  4  and  6  both  pressures  remained  at 
Pa/C  =  1-5  bar  and  the  relative  humidity  values  at  cathode  and  anode 
are  0.32,  0.40,  0.80,  0.92,  and  full  saturated  conditions  with  a  rel¬ 
ative  humidity  rha/c  =  1.0.  In  Tables  3  and  4  more  details  about  the 
experimental  conditions  are  given. 


Fig.  2.  Flow  chart  of  sequential  search  algorithm  for  deterministic  operating  area. 
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Fig.  3.  Flow  chart  of  sequential  search  algorithm  for  a  restricted  nondeterministic  operating  area. 
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Table  5 

Estimated  parameters  of  the  SEES  PEMFC  model  [8] 


Parameter 

Value 

Parameter 

Value 

li 

0.85620258 

& 

0.02675283 

h 

-0.05591375 

£9 

0.15494972 

h 

0.0005173 

£l0 

0.71394995 

$4 

0.04775405 

111 

0.05764422 

£5 

0.00430726 

Il2 

1.39237089 

Is 

12.5663035 

Il3 

4.48394466 

£7 

-0.02481325 

Il4 

0.00126953 

More  details  about  experimental  work  and  results  have  been 
previously  published  in  Refs.  [6,7]. 

3.2.  Fuel  cell  model 

In  a  previous  publication  of  our  team  [8]  a  valid  semi-empirical 
expert  system  (SEES)  fuel  cell  model  was  presented.  The  stationary 
cell  voltage  is  described  by  the  following  equation: 

Uce\\  =  Uocv  —  E^act  —  Ef0hm  —  Ud  if  —  E^hum-  (8) 

The  open  circuit  voltage  U0Cv  in  Eq.  (8)  is  described  as  a  function  of 
hydrogen  and  air  pressure.  Both  influences  are  lumped  together  in 
a  simple  equation  (Eq.  (9))  with  just  three  free  parameters: 

Uocv  =  §1  +  §2  ln(Pc)  +  §3  ln(Pa).  (9) 

If  the  cathode  pressure  rises,  the  first  term  in  Eq.  (10)  will  decrease 
leading  to  lower  activation  loss  Uac t>  The  logarithmic  dependency 
of  the  load  i  delivers  even  for  a  load  of  zero  valid  values  because  of 
a  leakage  current  parameter  §6. 

EEict  =  [§4-£5  ln(Pc)]  ln(i  +  §6)-  (10) 

The  quasi-linear  ohmic  resistance  part  of  the  current-voltage  curve 
shows  a  dependency  on  cathode  pressure  which  is  described  with 
a  linear  interaction  term  of  current  density  and  cathode  pressure 


as  given  in  the  following  equation: 

E Johm  =?7  i  +  ?8  i-Pc  (11) 

To  model  the  diffusion  loss  Udif,  Eq.  (12)  is  used. 

tfdif  =  §9Pc_1  exp(£10i).  (12) 

Eq.  (13)  calculates  the  last  term  of  Eq.  (8),  the  loss  due  to  insufficient 
membrane  humidity: 

Uhum  =  §11  exp(|12i)(1.01  -  rha)tl3  +§i4rha.  (13) 

The  given  PEMFC  model  uses  fcmax  =  4  inputs  (load  t,  cathode  pres¬ 
sure  pc,  anode  pressure pa,  and  anode  gas  relative  humidity  rha)  and 
0  =  14  model  parameters.  A  valid  set  of  these  parameters  is  given  in 
Table  5. 

Using  the  PEMFC  model  [Eqs.  (8)— (13)]  in  Eqs.  (2)  and  (3)  the 
FIM  can  be  calculated  and  used  for  the  introduced  design  criteria 
[Eqs.  (4)-(7)]. 

4.  Optimization  algorithm 

4. 1.  Algorithm  for  deterministic  operating  areas 

From  the  mathematical  point  of  view,  a  nonlinear  optimization 
problem  has  to  be  solved.  However,  the  design  task  is  restricted 
regarding  limitations  of  experimental  equipment  and  examinee. 
Another  challenge  is  the  large  number  of  degrees  of  freedom  that 
can  be  adjusted  to  optimize  the  design  criterion.  Looking  at  the  used 
fuel  cell  model  with  /<max  =  4  inputs  (0)  and  a  wanted  DoE  with 
0  =  50  independent  test  points,  200  parameters  would  influence 
the  design  criterion.  All  these  parameters  have  to  be  investigated 
to  get  the  optimal  experimental  design.  Most  used  numerical  gra¬ 
dient  methods  or  gradient  free  optimization  algorithms  will  have 
problems  with  disadvantageous  start  values  [9].  However,  genetic 
algorithms  are  less  sensitive  due  to  start  values  but  convergence 
speed  is  low  [10].  Therefore  the  optimal  DoE  is  a  promising  way. 

Since  the  A-,  E-  and  M-design  criteria  are  just  approximations 
of  the  searched  volumes  the  D-criterion  is  always  the  best  possible 


Table  6 

Candidate  set  for  the  simple  example 


# 

i  (A/cm2 ) 

Pa  (bar) 

# 

i  (A/cm2 ) 

Pa  (bar) 

# 

i  (A/cm2 ) 

Pa  (bar) 

# 

i  (A/cm2 ) 

Pa  (bar) 

1 

0 

1.1 

29 

0 

1.5 

57 

0 

1.9 

85 

0 

2.3 

2 

0.01 

1.1 

30 

0.01 

1.5 

58 

0.01 

1.9 

86 

0.01 

2.3 

3 

0.05 

1.1 

31 

0.05 

1.5 

59 

0.05 

1.9 

87 

0.05 

2.3 

4 

0.1 

1.1 

32 

0.1 

1.5 

60 

0.1 

1.9 

88 

0.1 

2.3 

5 

0.2 

1.1 

33 

0.2 

1.5 

61 

0.2 

1.9 

89 

0.2 

2.3 

6 

0.4 

1.1 

34 

0.4 

1.5 

62 

0.4 

1.9 

90 

0.4 

2.3 

7 

0.6 

1.1 

35 

0.6 

1.5 

63 

0.6 

1.9 

91 

0.6 

2.3 

8 

0.8 

1.1 

36 

0.8 

1.5 

64 

0.8 

1.9 

92 

0.8 

2.3 

9 

1 

1.1 

37 

1 

1.5 

65 

1 

1.9 

93 

1 

2.3 

10 

1.2 

1.1 

38 

1.2 

1.5 

66 

1.2 

1.9 

94 

1.2 

2.3 

11 

1.4 

1.1 

39 

1.4 

1.5 

67 

1.4 

1.9 

95 

1.4 

2.3 

12 

1.6 

1.1 

40 

1.6 

1.5 

68 

1.6 

1.9 

96 

1.6 

2.3 

13 

1.8 

1.1 

41 

1.8 

1.5 

69 

1.8 

1.9 

97 

1.8 

2.3 

14 

2 

1.1 

42 

2 

1.5 

70 

2 

1.9 

98 

2 

2.3 

15 

0 

1.3 

43 

0 

1.5 

71 

0 

2.1 

99 

0 

2.5 

16 

0.01 

1.3 

44 

0.01 

1.7 

72 

0.01 

2.1 

100 

0.01 

2.5 

17 

0.05 

1.3 

45 

0.05 

1.7 

73 

0.05 

2.1 

101 

0.05 

2.5 

18 

0.1 

1.3 

46 

0.1 

1.7 

74 

0.1 

2.1 

102 

0.1 

2.5 

19 

0.2 

1.3 

47 

0.2 

1.7 

75 

0.2 

2.1 

103 

0.2 

2.5 

20 

0.4 

1.3 

48 

0.4 

1.7 

76 

0.4 

2.1 

104 

0.4 

2.5 

21 

0.6 

1.3 

49 

0.6 

1.7 

77 

0.6 

2.1 

105 

0.6 

2.5 

22 

0.8 

1.3 

50 

0.8 

1.7 

78 

0.8 

2.1 

106 

0.8 

2.5 

23 

1 

1.3 

51 

1 

1.7 

79 

1 

2.1 

107 

1 

2.5 

24 

1.2 

1.3 

52 

1.2 

1.7 

80 

1.2 

2.1 

108 

1.2 

2.5 

25 

1.4 

1.3 

53 

1.4 

1.7 

81 

1.4 

2.1 

109 

1.4 

2.5 
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Fig.  4.  SSE  of  the  model  in  dependency  of  the  determinant  of  the  covariance  matrix 
of  the  used  experimental  design. 


choice.  The  other  criteria  are  just  advisable  for  very  fast  optimiza¬ 
tions  because  of  a  reduced  computing  time. 

If  there  is  only  a  low  hysteresis  in  fuel  cell  behaviour,  the  order  of 
performing  the  experiments  will  not  influence  much  the  measured 
voltage  and  it  is  possible  to  change  the  order  of  the  experimental 
conditions.  That  offers  the  opportunity  to  optimize  the  experi¬ 
mental  design  sequentially.  In  Fig.  2  the  flow  chart  of  a  developed 
sequential  search  (SS)  algorithm  which  identifies  the  best  0  exper¬ 
iments  of  MCImax  given  candidates  is  shown.  At  first,  the  number  of 
general  iterations  GImax  and  a  start  experimental  design,  randomly 
chosen  or  predetermined,  must  be  defined.  Then  the  first  row  (d  =  1 ) 
of  the  start  experimental  design  is  substituted  by  the  first  candidate. 
After  this  step  the  FIMi  is  determined  and  the  determinant  D\}  as 
in  the  following  equation: 

D,  =JD^  =det(FIM1/-FIM1).  (14) 

is  calculated.  Whereas  Jdmci  (Eq.  (5))  needs  to  be  minimized  (with 
MCI  as  the  run  index  in  the  range  of  1  to  MCImax),  DMci  must  be  max¬ 
imized  to  find  an  improved  experimental  design.  The  advantage  of 
using  DMci  instead  of  JoMa  is  the  avoidance  of  numerical  problems, 
caused  by  the  inverting  of  a  singular  matrix,  if  a  disadvantageous 
experimental  design  is  tested.  Both  steps  are  repeated  for  all  candi¬ 
dates.  At  the  end,  the  candidate  that  delivered  the  maximum  DMC i 
is  stored  and  taken  at  the  new  first  row.  Now  the  algorithm  contin¬ 
ues  with  the  second  row  (d  =  2)  and  again  all  candidates  are  tested 


Fig.  5.  SSE  of  the  model  in  dependency  of  the  ratio  of  experiments  (0)  to  the  number 
of  parameters  ((f)). 


Table  7 

Results  for  maximization  of  D  by  using  different  optimization  algorithms  at  a  work¬ 
station  (DELL  Precision  390  Workstation,  Intel  core  2  duo  1.86 GHz,  8  GB  RAM, 
WinXP64,  MatlabR2007b) 


SS 

GA 

LM 

Duration 

175  s 

9.5  h 

1.3  h 

D 

1.62  x  IQ”6 

6.53  x  IQ”8 

3.98  x  IQ”14 

and  the  candidate  with  the  largest  DMC i  is  stored.  This  is  repeated 
until  all  d  experiments  are  optimized  once  ( d  =  0 ).  After  that,  this 
loop  is  repeated  GI  max-times. 

4.2.  Algorithm  for  restricted  nondeterministic  operating  areas 

In  many  cases  no  experimental  data  will  be  available.  The  advan¬ 
tage  is  an  improved  quality  of  the  parameter  estimation  since  all 
operating  points  in  the  restricted  area  can  be  investigated  in  con¬ 
sideration  of  maximizing  the  covariance  matrix.  The  disadvantage 
is  that  our  sequential  search  has  to  be  upgraded  by  the  generation  of 
MCImax  candidates.  Fig.  3  shows  the  varied  flow  chart.  The  only  dif¬ 
ference  to  the  hitherto  used  algorithm  is  that  fcmax  random  numbers 
in  the  range  of  0-1  are  generated.  The  current  density  (9dl  is  gener¬ 
ated  randomly  by  a  Monte-Carlo  method  and  the  other  inputs  Gdk 
(/<  =  2. .  ./<max)  are  interpolated  in  the  nonlinear  restricted  operating 
area  in  dependence  of  the  current  density. 

The  generated  candidates  are  tested  and  compared  in  view  of 
the  parameter  estimation  quality  and  the  experimental  design  is 
gradual  improved  equal  to  the  previous  algorithm. 

5.  Simple  example  using  sequential  search  to  create 
D-optimal  DoE 

In  this  section  a  simple  example  is  given  demonstrate  how  to 
create  a  D-optimal  DoE  for  a  PEMFC  using  the  sequential  search 
approach. 

To  keep  the  example  clearly  represented,  only  the  influence  of 
anode  pressure  dependency  is  considered  in  this  section.  The  sim¬ 
plified  model  is  set-up  by  reducing  Eqs.  (8)-(12)  to  current  density 
i  an  anode  pressure  pa  dependent  terms.  The  resulting  equation 
reads  as 

Dceii  =  oi\  +  0:2  ln(pa)  -  c^3  ln(i  +  a 4)  -  a^i  -  a 6  exp^z).  (15) 

This  simple  model  take  /<max  =  2  inputs  into  account.  The  first  input 
( 6>i )  is  the  current  density  z,  and  the  second  ( ®u )  is  anode  pressure. 
Further  the  model  includes  in  total  0  =  7  model  parameters.  Thus, 
the  parameter  index  0  counts  from  1  to  7. 

Next  the  general  fimed  is  calculated.  For  an  arbitrary  set  of  exper¬ 
imental  conditions  d,  fim  includes  7  columns  and  reads  as  given  in 
the  following  equation: 

"1 

In  (p„) 

-ln(i  +  a4) 


-z 

-  exp  (a7  •  z) 

_-oi6  ■  i  •  exp  (a?7  •  z)_ 

To  calculate  the  design  criterion  D,  as  is  given  in  Eq.  (14),  the  initial 
model  parameter  must  be  given  to  Eq.  (16).  Each  a  corresponds 
to  a  §  of  the  original  model,  as  described  in  Eqs.  (8)— (11 ),  which 
are  given  in  Table  5,  and  are  calculated  with  an  assumed  cathode 


-1  T 


(16) 
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Table  8 

Calculated  determinant  by  sequential  search  after  each  iteration  step  for  different  start  sets 


GI 

D 

0 

1 

2 

3 

4 

5 

Start  set  A 

0 

1.45  x  10“6 

1.61  x  10“6 

1.62  x  10“6 

1.62  x  10“6 

1.62  x  10“6 

Start  set  B 

2.62  x  10-18 

1.58  x  10“6 

1.61  x  10“6 

1.62  xlO-6 

1.62  xlO-6 

1.62  xlO”6 

Start  set  C 

2.06  x  10“16 

1.57  x  10“6 

1.61  x  10“6 

1.62  x  10“6 

1.62  x  10“6 

1.62  x  10-6 

pressure  of  1.5  bar  following  equation: 


r§1+§2.ln(1.5)l 

?3 

§4  -  ?5  •  ln(1.5) 
§6 

?7+?8-1.5 

^g.1.5-1 

_£l0 


(17) 


Now  we  assume  Ui-curves  are  experimentally  measured  from  0  up 
to  2  A/cm2  with  14  load  points  each  Ui-curve,  and  anode  pressure 
is  varied  in  the  range  from  1.1  to  2.5  bar  in  0.2  bar  steps. 

To  design  the  D-optimal  DoE  the  maximal  number  of  exper¬ 
iments  0  has  to  be  defined.  In  this  example,  two  times  more 
experiments  shall  be  carried  out  to  identify  model  parameters  well. 
This  leads  to  0  =  14.  A  traditional  DoE  set-up  by  a  scientist  would 
now  measure  for  every  anode  pressure  a  single  Ui-curves,  which 
would  lead  to  8  x  14  =  112  single  measurements. 

These  112  experimental  conditions  are  given  in  Table  6,  and 
are  used  as  candidate  set  (MCImax  =  112)  for  the  sequential  search 
algorithm. 

The  sequential  search  algorithm  is  initialized  with  the  first 
14  candidates,  which  equals  a  initial  design  criterion  value  of 
Dinitiai  = -1.77  x  10-44.  This  very  poor  value  is  related  to  the  fact 
that  there  are  only  initial  candidates  atpa  =  1.1  bar.  During  first  run 
(d  =  1 )  of  the  algorithm  the  first  entry  of  initial  DoE  is  sequentially 
substituted  by  all  112  candidates,  and  the  resulting  design  criteria 
are  calculated.  At  the  end  candidate  #108  gives  best  design  crite¬ 
rion  for  first  entry  with  a  D^  =  1.46  x  10-28,  and  therefore,  is  stored 
to  the  DoE.  The  design  criterion  improved  by  72  orders  of  magnitude 
just  by  replacing  the  first  entry.  Best  candidate  for  second  entry  is 
found  to  be  candidate  #100,  with  D2  =  2.29  x  10-28.  After  complete 
optimization  of  all  14  entries  in  the  DoE  the  sequential  search  algo¬ 
rithm  repeats  two  times  (GImax  =  3)  the  complete  process.  The  final 
design  criterion  (Dfinal)  can  be  improved  to  Dfinal  =  1.84  x  10-27,  and 
contains  the  candidates  #1,  5,  8,  10,  13,  14,  99,  100,  103,  105,  106, 
109,  111,  and  122. 


6.  Experimental  validation 

6.1.  D-criterion 


For  the  validation  of  the  D-optimal  design  criterion  the  0  =  14 
parameters  of  the  PEMFC  model  are  estimated  for  different  sets  of 
0  =  50  data  points  selected  from  the  405  experimental  data  points. 
For  the  evaluation  of  the  estimated  model,  the  sum  of  squares  due 
to  errors  (SSE)  between  measured  Ud  and  model  predicted  fuel  cell 
voltage  Ud,  as  given  in  the  following  equation,  is  used. 

0 

SSE  =  ^(Ud-0d)2  (18) 

d= 1 

The  SSE  of  each  validation  trail  (SSEvalidation)  is  plotted  in  Fig.  4 
against  the  achieved  determinant  (Dvalidation)  to  observe  potential 
correlations.  A  suitable  scaling  of  the  respective  optimum  (SSEmin 
and  Dmax)  is  chosen  for  a  better  illustration. 


Clearly  a  strong  correlation  between  both  states  can  be  seen. 
Moreover  a  strong  decrease  of  SSE  can  be  seen  with  increased  deter¬ 
minant.  A  model  with  parameters  from  an  experimental  design 
with  a  small  determinant  has  a  significantly  bigger  SSE  than  a 
model  with  parameters  from  an  experimental  design  with  a  large 
determinant.  Since  the  SSE  decreases  strongly  exponential  with 
an  increasing  determinant,  the  D-criterion  is  for  sure  an  excellent 
method  to  evaluate  the  quality  of  parameter  estimation  and  was 
used  for  all  further  validations. 

6.2.  Sequential  search 

The  performances  of  three  different  optimization  approaches  to 
maximize  D  were  investigated.  The  proposed  sequential  search  was 
compared  to  a  Levenberg-Marquadt  (LM)  algorithm  and  a  genetic 
algorithm  (GA).  Further  information  about  both  methods  can  be 
found  in  Refs.  [9,11,12].  All  three  algorithms  were  initialized  with 
the  same  set  of  randomly  chosen  50  experimental  conditions  with 
Djnitiai  =  2.06-16.  In  Table  7  the  maximum  achieved  determinant  and 
the  duration  to  get  this  result  is  given  for  each  algorithm. 

Investigating  various  configurations  of  the  LM  algorithm,  only 
a  slight  improvement  of  D  could  be  achieved  and  1  h  20  min  were 
needed  to  calculate  the  best  result.  In  most  of  the  investigated  LM 
algorithm  set-ups  even  no  improvement  of  the  design  criterion  was 
achieved.  Using  a  GA  an  improved  maximization  by  8  orders  of 
magnitude  of  D  is  gained  by  cost  of  a  very  long  calculating  time. 
But  also  the  GA  depends  strongly  on  the  used  algorithm  options 
and  delivered  no  results  in  many  cases.  Using  the  sequential  search 
algorithm  a  24.8  times  greater  Dina  very  fast  calculation  of  less  than 
3  min  was  found  at  a  DELL  Precision  390  Workstation,  using  Intel 
core  2  duo  1.86  GHz  P,  8  GB  RAM,  WinXP64,  and  MatlabR2007b. 

6.3.  Robustness  of  sequential  search 

Usually  there  are  no  good  start  values  available,  which  mean  a 
good  experimental  design.  Different  sets  were  used  for  the  sequen¬ 
tial  search  of  the  optimal  50  candidates  for  investigating  the  quality 
of  this  method  in  consideration  of  the  start  values.  The  results  for 
three  different  start  sets  are  displayed  in  Table  8.  Even  starting  with 
a  very  poor  experimental  design  as  start  set  A,  the  algorithm  works 
properly.  The  maximum  determinant  was  found  for  every  start  set 
already  after  three  iterations  (GI  =  3)  and  even  the  first  and  second 
iteration  delivered  useful  determinants.  Thus  the  sequential  search 
is  not  only  a  very  fast,  but  also  a  very  robust  method  to  optimize 
the  proposed  design  criterion. 

6.4.  Experimental  effort 

Another  issue  of  this  work  is  the  question  how  many  exper¬ 
iments  are  needed  for  getting  satisfying  parameter  estimations. 
Therefore  the  ratio  of  the  number  experiments  to  the  number 
parameters  is  plotted  against  the  SSE  between  model  and  exper¬ 
imental  data  in  Fig.  5.  A  significant  improvement  of  the  model  can 
be  observed  till  a  ratio  of  four.  More  experiments  deliver  just  a  low 
improvement  and  the  SSE  reaches  a  steady  state  which  can  not  be 
decreased  because  of  inaccuracies  of  the  model  itself. 
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7.  Discussion 

A  constraint  in  this  work  is  that  no  dynamic  database  was  inves¬ 
tigated.  Furthermore,  the  gained  experimental  design  can  be  a  local 
optimum  and  so  there  is  no  guarantee  that  not  a  higher,  unnoticed 
global  optimum  exists. 

For  reducing  the  computation  time  the  replacement  of  the  deter¬ 
mined  number  of  iterations  a  in  Fig.  2  by  an  additional  algorithm 
stop  similar  to  the  abort  criterion  for  numerical  algorithm  is  one 
of  the  next  steps.  If  gradient  of  D  does  not  improve  significantly 
between  two  loops,  no  further  iterations  are  needed.  The  advan¬ 
tage  would  be  an  optimized  number  of  iterations  and  the  avoiding 
of  unnecessary  steps  as  seen  in  Table  8. 

Also  interesting  would  be  a  hybrid  approach  for  optimizing  the 
experimental  plan.  The  sequential  insertion  of  all  candidates  or 
the  generation  of  random  variables  could  be  replaced  by  a  local 
search  of  the  optimum  row  by  using  a  numerical  algorithm  like  for 
example  a  LM  algorithm. 

Problems  with  the  introduced  approach  could  appear,  when  old 
stacks  with  thick  membranes  are  used.  These  stacks  often  have  a 
hysteresis  of  about  20  mV.  In  this  case,  the  order  of  experimental 
conditions  is  not  negligible  and  the  SS  may  deliver  invalid  results. 

8.  Conclusions 

It  is  shown  that  the  model  is  valid  for  DoE  and  prior  improve¬ 
ment  of  the  experimental  design  is  a  promising  way  for  enhanced 
parameter  estimation.  Especially,  the  D-optimal  design  approach  is 
suitable  to  identify  model  parameters  efficiently  as  Fig.  4  demon¬ 
strates.  The  developed  sequential  search  is  compared  to  usual 
used  methods  of  optimization  and  shows  a  significant  lower  cal¬ 
culating  time  as  well  as  a  couple  of  orders  of  magnitude  larger 
determinant. 

The  sequential  search  is  moreover  recommendable  for  DoE 
because  of  a  low  number  of  iterations  are  needed  and  the  method 
proves  a  high  robustness  against  disadvantageous  chosen  start 
designs.  Already  a  number  of  experimental  data  points,  which  is 
as  big  as  four  times  the  number  of  model  parameters,  delivers  an 
accuracy  model  simulation  result,  further  experiments  delivers  just 
low  improvements  and  are  in  the  most  cases  not  recommendable 
in  view  of  increasing  experimental  effort. 

All  in  all  the  prior  optimization  of  the  experimental  design  in 
order  to  improve  the  parameter  estimation  is  an  intelligent  way 
to  save  time  and  thus  to  reduce  the  cost  of  experimental  work 
significantly. 

The  benefit  performing  such  a  DoE  is,  to  use  expensive  exper¬ 
imental  resources  much  more  efficiently.  Traditionally,  fuel  cells 


are  characterized  by  Ui-curves,  as  described  in  Refs.  [6,7].  There, 
the  characterization  for  just  4  inputs  causes  an  effort  of  405  sin¬ 
gle  measurements.  Using  the  proposed  D-optimal  DoE  together, 
with  the  sequential  search,  this  effort  can  be  reduced  by  a  factor 
of  more  than  8  to  a  total  cost  of  50  single  measurements.  With  the 
proposed  method,  experimental  groups  are  able  to  characterized 
fuel  cells  or  fuel  cell  stacks  in  a  very  large  operation  area  with  an 
affordable  experimental  effort.  This  helps  to  speed  up  global  fuel 
cell  modelling  tremendously,  and  provide  urgent  model-input  to 
fuel  cell  system  concept,  design,  and  engineering  work.  This  con¬ 
tribution  helps  to  provide  efficiently  high  quality  models,  which  are 
the  key-tools  in  modern  drive  train  development.  Reducing  engi¬ 
neering  costs,  and  improving  fuel  cell  drive  train  performance  is  the 
major  hurdle  in  commercialization  this  zero  emission  technology, 
and  is  supported  by  this  work. 

Further  work  is  advisable  in  terms  of  taking  over  this  general 
method  to  other  fields  of  science  and  so  to  establish  our  method  as 
a  general  valid  way  for  optimizing  experimental  designs. 
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